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ABSTRACT 


In this investigation, two-dimensional turbulent flow of incom- 
pressible Newtonian fluids in sharply curved 180° turn around ducts is 
studied. Results of an approximate numerical flow field analysis 
utilizing an orthogonal, body-fitted, curvilinear coordinate system are 
compared to results based on a traditional cylindrical reference frame. 
Qualitative indication of general streamfield characteristics as well 
as quantitative benchmarks for the planning of future experimentation 
are provided. In addition, preliminary results of an augmented k-e 
turbulence model analysis, which explicitly accounts for the effects of 
streamline curvature and pressure strain in internal turbulent flows, 
are presented. Specific model difficulties are discussed and compari- 
sons with standard k-e model predictions are included. 
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NOMENCLATURE 


Symbol 


Definition 



constant (=0.09) in traditional k-e model turbulent 
viscosity expression 

spatially variant function in modified k-E model turbulent 
viscosity expression (u t =C pk 2 /e) . 

constant m wall function correction (=</C 3 n ‘*), see 
equation (13) ^ 

duct width 

dissipation term in the Reynolds stress transport equations 

east side cell, see Figure 12 

cubic polynomial in v^CP, see equation (12) 

turbulent kinetic energy 
turbulence length scale 

turbulent length scale for flow over a plane 

north side cell, see Figure 12 
time averaged pressure 
turbulent energy production rate 

exit plane pressure assumed uniform 

production term in the Reynolds stress transport equations 

radial direction coordinate 

duct inner convex wall radius of curvature 

duct outer concave wall radius of curvature 

duct centerline radius of curvature (=[r.+r 1/2) 

l o J 

south side cell, see Figure 12 

momentum source terms, see equations (5) and (6) 


u r Ug, T7, ^ Reynolds stress components in cylindrical coordinates 
U time averaged streamwise velocity 


u e 


w 

y 


cross-stream averaged streamwise velocity 

time averaged velocity components in cylindrical 

coordinates 

west side cell, see Figure 12 

distance from the duct inner convex surface measured in the 
radial direction 
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Definition 


Symbol 

8 empirical constant in curvature modified length scale, 

see equation (1) 

£ energy dissipation rate 

dissipation term in the Reynolds stress transport equations 
K von Karman constant 

U absolute viscosity 

Ut turbulent viscosity 

heff effective viscosity (=u t +u) 

Ilij pressure strain term in the Reynolds stress transport 

equations 

p mass density 

0 angular position measured from the duct entrance plane 

BFC body-fitted coordinates 

FDP fully developed pipe or duct 

SSME space shuttle main engine 

TAD turn around duct 
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INTRODUCTION 


Performance demands of modern machinery require designs exhibiting 
both dimensional and operational economy. In applications which require 
the motion of a working fluid, these economic considerations have fueled 
a proliferation of devices containing complex flow passages through which 
fluids travel at extreme rates. Examples include flows in turbomachinery, 
nozzles, diffusers, and curved ducts. Of primary interest in this inves- 
tigation are flows in sharply curved turn around ducts (TAD's) common in 
rocket engine design, including current and projected versions of the 
space shuttle main engine (SSME). 


Despite the design impetus provided by increased performance 
requirements, a detailed understanding of turbulent fluid motion with 
significant streamline curvature is inhibited by formidable obstacles to 
both experimental measurement and computational study. Since the system 
of motion governing equations for turbulent flows is essentially incom- 
plete, advances in the computational prediction of turbulent field char- 
acteristics are tied to turbulence models requiring a sound empirical 
basis. In the case of sharply curved passages and extreme flow rates 
prevalent in modern equipment, the experimental basis required for 
detailed flow prediction is difficult to obtain with confidence. 


An example of a high rate, strongly curved flow field with little 
experimental basis on which to develop adequate turbulence models occurs 
in the SSME fuel side turbopump TAD. A simplified two-dimensional 
version of this TAD configuration is depicted in Figure 1. Although 
several studies of turbulent curved duct flow have been reported [1-7], 
none deals with curvatures or rates of the magnitude occurring in the 
SSME. The state of turbulence modeling in general is such that little 
confidence can be placed in extrapolation of the results of these inves- 
tigations to a broader class of curved flows. 


The sensitivity of turbulent flow characteristics to even mild 
streamline curvatures has been discussed by Bradshaw [8]. In this 
extensive review, a curvature modified length scale of the form 


Z 


(1±S 


U 0 /r 



)lo 


( 1 ) 


was proposed. In this relation 8 is an empirical constant of order 10, 

l 0 is a length scale appropriate for plane flow, and subscripts following 

a comma indicate differentiation with respect to the specified variable 

(e.g., U Q = 3U Q /3r). Inclusion of this curvature modified scale into a 
° 0 r 0 
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general mixing length model of turbulent motion is hampered by the flow 
dependence of both 3 and % 0 . Lack of extensive applicability is charac- 
teristic of all mixing length modeling strategies. 

Launder et al. [4] have proposed a two equation closure model 
with the length scale dependence obtainecU from solution of the equations 
for turbulent kinetic energy and the energy dissipation rate. The 
direct effect of curvature is introduced as a single term in the dissi- 
pation equation which is proportional to the turbulent Richardson number. 
The required proportionality constant appears as a new empiricism. 

Recently, an extended k-e turbulence model, which uses algebraic 
approximations for the Reynolds stress transport equations in order to 
develop an explicit curvature dependent expression for in = CMc 2 /e, 

has been reported by Pourahraadi and Humphrey [9,10] and applied to 
curved duct flows. The new expression for Cy accounts simultaneously 
for effects of streamline curvature and pressure strain with wall 
dampening. The development is a formidable task, calling upon a variety 
of traditional turbulence modeling approximations [11-15]. The approach 
taken by Pourahmadi and Humphrey presents a promising explicit curvature 
model with roots, albeit approximate, in the Reynolds stress transport 
equations. Unfortunately, previous results based on this computational 
procedure appear to be based on a rather fundamental algebraic error 
which will be addressed in subsequent sections of this report. 

Another approach which has been pursued is to include curvature 
effects in a modified wall law [3,16]. This method is well suited for 
near wall treatment in coarse grid numerical calculations of turbulent 
boundary layers under adverse pressure gradients. Incorporation into 
computational procedures involving strongly curved, fully elliptic flow 
fields is certainly plausible. 

Because of the scarcity of reliable experimental data involving 
high rate flows in sharply curved passages, and the difficulty in extrap- 
olating turbulent model predictions , none of the above described curva- 
ture modified turbulence models can be objectively recommended. Each, 
however, presents an opportunity for computational investigation, refine- 
ment, and comparison with results of standard turbulence models as well 
as experimentation. 
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OBJECTIVES 


The objectives of this investigation have a dual character. On 
one hand, computational predictions of turbulent flow properties in 
sharply curved TAD's are sought in order to provide additional infor- 
mation on the operational characteristics of devices incorporating such 
flows. Of specific interest are TAD flows in the SSME fuel side turbo- 
pump. As such, results of this investigation fall into the category of 
computational analysis. A second view is more appropriate, however, 
with results of the investigation taken as presenting contributions to 
the base of computational experimentation needed to improve general 
models of turbulent flow with significant streamline curvature. The 
data base building aspect of this study is considered to be the primary 
goal. 


In order to achieve this primary goal, three specific objectives 
were established. 


1. To estimate the general characteristics of turbulent flows 
in sharply curved TAD's using a standard k-e turbulence 
model on an abbreviated geometry as displayed in Figure 2. 

2. To implement and test a specific turbulence model which 
explicitly incorporates the effects of streamline curvature 
and pressure strain. 


3. To examine the elliptic character of the governing motion 
equations by implementing an orthogonal, body-fitted 
coordinate reference frame analysis on an extended duct 
geometry shown in Figure 3. 


Background information concerning each of the three stated objectives 
is provided in the following section of this report. 
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Figure 3. Geometry and grid pattern for BFC TAD analyses: 20 X 80 grid. 


BACKGROUND INFORMATION 


Equations governing the motion of*steady, two-dimensional, tur- 
bulent flow of an incompressible, Newtonian fluid in cylindrical coor- 
dinates are given below [10]. 


Continuity : 
r-Momentum: 


0-Momentum: 


U r,r + u e,e /r + V - 0 


p[U r,r + Vr,e /r ' U 6 /r ' 


( 2 ) 


- P ,r + (ll eff rU r,r ) ,r /r + ' u aff U r,e /r) ,0 /r 


' u eff U r /r! ' 2 u eff U e,e /ra + S r 
p[U r U 9i r + U e U ee /r ♦ UUg/rJ 

" - P ,8 /r + <M eff r Vr ) ,r /r 


(3) 


+ (M eff U 6,6 /r) ,e /r “eff'V 1 " 

+ 2lJ eff u r,e /r! + S 6 


( 4 ) 


In equations (3) and (4) above, the Boussmesq assumption has 
been employed to model the Reynolds stresses in terms of mean velocity 
gradients and a turbulent viscosity component. The momentum source 
terms S r and Sq are given by the relations below. 


S = 


= 


>V ( V r) ,rJ,e /r + <V U r,r > ,r /r ‘ W' 1 * 

[U t (2U./r <- U 6f9 /r)] i9 /r + (u t (U r 0 - V /r J, r 

+ U t (U e,r " V r)/r 


(5) 


( 6 ) 


In order to provide and solve for spatial variation of the tur- 
bulent viscosity U t , the k-e turbulence model as described by Launder 
and Spalding [17] was employed. This two equation model has been com- 
pared favorably to a variety of turbulence modeling schemes [18]. It 
has been tested on a wide range of flow configurations, exhibiting some 
degree of universality with an acceptable degree of complexity. Fol- 
lowing this approach, is determined by the relation 


= c y pkZ/E ’ 


(7) 
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where the turbulent kinetic energy and dissipation rate are determined 
upon solving the differential transport equations given below [10,19]. 


turbulent kinetic energy k: 

p[U r\r + U e k ,e /r] * ^ p eff sk ,r /o k^ ,r /r 


+ (u eff k i e /o ic ) I 0 /r2 “ p£ + P d 


( 8 ) 


energy dissipation rate e: 


0tU r E ,r + V,0 /r > ' <“eff rE ,r /o e >,r /r + <U ef f E ,e> 


+ C £l eP d /k - C £2 pe 2 /k 


(9) 


The production term in equations (8) and (9) is given by the relation 


P d - u t (2[(« rjr ) a ♦ (U e y> ! - V u r,e /r + V )/r 

+ W + 2U e,e> /ra + u r,e u e,r /r l + ( V r)a + <Vr> 2 
+ (U r yr) 2 } . (10) 

The numerical procedure used to solve the discretized form of 
equations (2), (3), (4), (8) and (9) was provided by the proprietary 
PHOENICS computer code supplied and supported by CHAM of North America. 
A discussion of the philosophical and mathematical basis of this code 
is provided in reference [19]. A good presentation of the numerical 
procedures associated with a large family of CFD codes including the 
PHOENICS code is available in the text by Patanker [20], 

The standard k-e turbulence model, the basis of which is pre- 
sented above for cylindrical coordinate systems, was employed in compu- 
tational analyses of simple geometry (Figure 2) TAD flows. A 20 x 30 
grid was employed in these analyses as displayed between planes BB and 
CC of Figure 3. In all simple TAD configuration analyses, inlet plane 
velocity, turbulent kinetic energy, and dissipation rate conditions 
were specified. A constant exit plane (CC in Figure 3) pressure was 
prescribed and boundary condition closure was achieved by assuming zero 
velocity field gradients in the exit plane. 

The Pourahmadi and Humphrey [9,10] extension of the standard 
model was also implemented. As discussed above, this modified turbu- 
lence model incorporates the effects of streamline curvature and 
pressure strain in the fluid by using a variety of approximations to 
reduce the Reynolds stress transport equations to a set of algebraic 
equations. Because of the complexity of the algebraic manipulations 
employed, only a brief outline of the model derivation and implemen- 
tation is presented below. The reader is referred to the original work 
in reference [9] for a detailed development. 


XXXV-8 



The derivation starts from a high Reynolds number form of the 

u^u,- transport equation [9] characterized by the general form shown 
below. 


D Uj[Uj 
Dt 


= P. . + £ . + fl. . + D. . 
ij ij *ij ij 


( 11 ) 


The development then proceeds in the following stepwise manner. 


1) The u x Uj transport equations are reduced to a system 
of algebraic equations by invoking a variety of 
approximations including 

A) Rotta's [12] return to isotropy principle in Jlji 

B) Launder, Reece and Rodi's [13] expression for the 
contribution to II X ^ from interaction between mean 
strain and fluctuating velocities 

C) Daly and Harlow’s [14] correction to due to 
wall effects 

D) Gibson and Launder’s [15] correction to II due to 

wall effects 

E) Rodi's [11] assumption that u^u./k is constant 
throughout the flow field. 

2) The resulting system of algebraic equations is solved 
for the turbulent stresses. 

3) The Boussinesq approximation is then employed, trans- 
forming the solution into an algebraic relation for C^. 

4) The turbulent energy production term is rearranged into 
the form 

W e >Vr ’ 8<V £ ' V r ’ 

where g is a complex function of the quantities indicated, 
including the various partials of the mean velocity, j. 

5) The relation for (k/e)Ug r from step 4 is substituted 
into the algebraic relation for from step 3. Simpli- 
fication of the resulting expression leads to the cubic 
polynomial relation 

‘ * a i c u + * 2 C l' ! + a 3 - 0 (12) 


where a^a^a^ = a^ ,a 2 ,a 3 (P d /e, U 0 /r, f, Ih 

6) The wall correction function f presented in reference 
[13] and appearing in step 5 is modified to the form 


f = kiiA r I . (y/°) m i 
" C w £ 1 y ~D^y 1 


( 13 ) 
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thereby permitting asymmetric wall corrections to the pressure 
strain. In this expression, the asymmetry measure m appears 
as a new empirical constant. 

7) Following Bradshaw [8], the expression for dissipation at near 
wall nodes is modified by a multiplicative factor of 
l/[l±8Ug/(Ug r r)] in order to introduce the influence of extra 

strain curvature on the magnitude of the near wall turbulence 
length scale. This expression also contains a new empirical 
constant, 6. 

Suggested values for constants m and 8 are 7.95 and the order of 10 res- 
pectively. 

Although the complex development terminating at equation (12) is 
quite impressive, an unfortunate algebraic error is made in establishing 
the roots of this equation. Instead of a single positive real root for 
as predicted in references [9] and [10], multiple real roots, not 
necessarily positive, can (and did) occur in practice. The implemented 
version of the Pourahmadi and Humphrey modification included a corrected 
method of solving for the curvature dependent values of Adding to 

the difficulties associated with this implementation was the need to 
prescribe a selection criterion for regions in which multiple roots 
occur. This problem was never resolved with complete success as will be 
discussed more fully in the next section of this report. 

In order to more fully understand the elliptic character of devel- 
oping TAD flows, analyses utilizing an orthogonal body-fitted coordinate 
(BFC) reference system on an expanded duct geometry were performed. Two 
such analyses were performed; one on the full geometry, 20x80 grid 
system depicted in Figure 3, and a second on a reduced geometry 20x60 
grid depicted between planes AA and CC of the same figure. In each case 
the 20 x 30 grid configuration of the simple TAD analyses was preserved 
in the region bounded by surfaces BB and CC in Figure 3. In this manner 
comparison of results using both simple TAD configurations and expanded 
BFC geometries was facilitated. 

Numerical computations employing the grid configuration of Figure 
3 were performed utilizing the BFC option of the PH0ENICS computer code 
described in reference [21]. Although potentially more accurate non- 
orthogonal grid specifications could have been implemented within the 
framework of the PH0ENICS code, direct comparison with results of simple 
cylindrical system analyses would have then become more difficult. An 
excellent comprehensive review of methods to numerically generate curvi- 
linear coordinate systems with coordinate lines coincident with boundary 
segments is given m reference [22], 
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RESULTS 


Results of this investigation can be separated into three cate- 
gories dependent upon coordinate reference frame and turbulence model 
employed. Analyses 1, 2, and 6 described in Table I were performed on 
simple TAD geometries using a cylindrical coordinate system and a stan- 
dard k-e turbulence model [17]. In analysis 3 the modified k-e model of 
Pourahmadi and Humphrey [9,10], with the correction discussed in the 
previous section, was employed. The simple TAD geometry displayed in 
Figure 2 was maintained. Analyses 4 and 5 employed orthogonal BFC refer- 
ence systems with a standard k-e turbulence model . The grid system de- 
picted in Figure 3 was used in analysis 4. The same grid pattern, res- 
tricted to the corresponding analysis geometry, was employed in all 
cases. In addition, the flow Reynolds number was fixed at 10^ and a 
curvature ratio R c /D = 1 was prescribed. These values obviously refer 
to a high rate flow with extremely sharp curvature. A complete summary 
of the computational analyses of this study is presented in Table I. 

In Figures 4, 5, and 6 are displayed developing profiles of veloc- 
ity, pressure, and turbulent kinetic energy respectively, for the case 
of a simple TAD geometry with inlet plug flow. In response to the strong 
cross stream pressure gradients evident throughout the flow domain, the 
velocity profiles are observed to become increasingly asymmetric. The 
most substantial pressure adjustments occur within the final 30° of duct 
travel, and in response there is a significant acceleration of flow near 
the outer wall with a corresponding deceleration near the inner, convex, 
duct wall. The net effect is a nearly linear velocity profile very near 
the exit plane reminiscent of a solid body rotation. Quantitatively, the 
near exit plane velocity varies from just under 25 % to somewhat over 
140% of the average flow velocity. In the entrance region of the duct, 
the plug flow characteristic is observed to persist virtually unchanged 
past the midplane (9 = 90°) of the duct. The large cross stream pressure 
gradient predicted for the TAD inlet region together with the specified 
plug flow entrance velocity profile is suggestive of a very large swirl 
chamber feed to the TAD. This type of configuration is of course not 
prevalent in rocket engine or turbomachinery applications. 

The turbulent kinetic energy profiles displayed in Figure 6 ex- 
hibit rapid decay in the plug flow core region as expected. An initial 
decay in the near wall turbulent kinetic energy is followed by a signif- 
icant rise in the accelerating concave wall region, and a continuing 
gradual decay in the decelerating convex wall region. Near the exit 
plane, the turbulent energy profile develops the knee observed [1,2] and 
predicted [10] in previous curved duct investigations conducted at 
smaller values of Re and R c /D. 
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Table I 


Analyses Summary 


Analysis 

No. 

Reference 

Frame 

Figure 2 
I/O Plane 

Grid 

Inlet U 
Profile 

Turbulence 

Model 

1 

Cyl. 

BB/CC 

20x30 

Plug 

k-e 

2 

Cyl. 

BB/CC 

20x30 

FDP 

k-e 

3 

Cyl. 

BB/CC 

20x30 

Plug 

Modified 

k-e 

4 

BFC 

AA/EE 

20 x 80 

Plug 

k-e 

5 

BFC 

AA/CC 

20 x 60 

Plug 

k-e 

6 

Cyl. 

BB/CC 

20x30 

Plane BB 
from 5 

k-e 

Analyses Constants: Re = 10^ 

Turbulence Inlet Conditions 

5 R /D 

c 

— 9 

= 1.0 


For Plug Inlet U Profile 


k/U = 0.005 
e/(U 3 /D) = (0.005) 3/ 2 

/0.01 

Exit Plane 

Pressure: 

Uniform 

across exit 

plane 
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Results of a simple TAD geometry analysis with a fully developed 
turbulent duct flow inlet condition are displayed in Figures 7, 8, and 9. 
The general shape characteristics and scale of the pressure profiles in 
Figure 8 are similar to the previous plug flow inlet case. The develop- 
ment sequence for the velocity profiles is also similar to the analysis 
1 case with little profile adjustment prior to the 90° midplane of the 
TAD, and significant velocity gradient increases and profile lineariza- 
tion in the rapid adjustment exit region. The salient feature displayed 
in Figure 7 is the appearance of a recirculation region near the inner 
wall of the duct. This region extends back somewhat over 6° from the 
exit plane and occupies approximately 4 % of the duct cross section at 
the exit plane. The occurrence of a negative velocity region in the 
exit plane presents certain difficulties involving computational error 
as discussed m reference [21]. Despite this problem, flow separation 
from the inner convex wall is predicted to occur between 0 = 171° and 
0 = 174°. This phenomena is not observed to occur in less sharply 
curved channels with smaller values of the flow Reynolds number [1,2,6]. 

The turbulent kinetic energy profiles depicted in Figure 9 differ 
substantially from the plug flow inlet case both in shape and develop- 
ment history. The initial profile shape is symmetric and of the same 
magnitude as displayed for the analysis 1 entrance profile. Because of 
the entrance cross stream velocity gradients present m the FDP entrance 
case, however, the turbulent kinetic energy is predicted to increase 
significantly in the flow direction away from the convex wall. This is 
in sharp contrast to the rapid post entrance decay of k exhibited in 
Figure 6 for the plug flow inlet case. Near the exit plane, cross 
stream velocity gradients and the energy dissipation rate away from the 
convex wall increase rapidly in the flow direction, and the increase in 
turbulent energy is arrested. In the vicinity of the convex wall there 
is a substantial adverse pressure gradient leading to deceleration and, 
as the exit plane is approached, flow separation. This is accompanied 
by a sharp decrease in the inner wall region values of k. 

Direct comparison of the downstream velocity and pressure pro- 
files predicted in analyses 1 and 2 are displayed in Figures 10 and 11. 
The similarities in each of the profile pairs is obvious. The near exit 
plane pressure adjustment is predicted to occur more suddenly in the 
case of plug flow inlet conditions. It is evident that the inner wall 
characteristics of the near exit plane flow are substantially affected 
by changes in the inflow boundary conditions. 

In the third analysis of the computational sequence, the modified 
k-e model of Pourahmadi and Humphrey [9,10] was implemented in the user 
definition section of the PHOENICS routine. Severe difficulties were 
initially encountered in attempts to obtain reasonable solutions. Fur- 
ther investigation revealed the equation (12) rooting error described 
previously. 

In the context of the PHOENICS implementation, the Pourahmadi and 
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Humphrey [8,10] modification for a specific finite difference cell P re- 
quires manipulation of the staggered grid velocity vectors indicated in 
Figure 12. In addition, grid distances from each of the indicated neigh- 
boring cells to the P cell boundaries are needed together with estimates 
of k and £ at P. The entire manipulation is computationally intense, 
often tripling the clock time of analysis.* 

Despite the increased computational effort, convergence diffi- 
culties were repeatedly encountered. The source of the problem was the 
behavior of the polynomial F(i/C^), the roots of which indicate candidate 
values of Cy. Typical F(/C^) loci are depicted in Figure 13. These 
polynomial traces were arrived at by computationally testing the modified 
k-e turbulence model on a simple TAD geometry, plug inlet flow problem. 
Analysis 3 as described in Table I was initiated at the standard k-£ 
model solution value of analysis 1. An intermediate flow field result 
was determined by arbitrarily selectin g th e principal theta transfor- 
mation root of the cubic polynomial F(/Cy) whenever multiple roots oc- 
curred. The trigonometric or theta transformation method for determining 
roots of a cubic polynomial is described in most mathmatical handbooks 
(see e.g. [23]). 

It is evident from Figure 13 that for 0 =168°, the near concave 
wall behavior of F was that of a cubic polynomial with multiple real 
roots. This behavior was not isolated, occurring irregularly throughout 
the flow field and frequently in near wall regions. Unfortunately it is 
in these near wall regions that the modified k-£ model was expected to 
yield its most meaningful results [9,10]. The behavior of the F(/C^) 
polynomial at y/D = 0.9875 (the node closest to the concave wall) was 
particularly disturbing. As indicated m Figure 13, the trace of F was 
virtually colinear with the Cy axis for -0.2<Cy<0.2. Three distinct real 
roots were determined at this nodal position, however, the ill natured 
character of F is apparent, making the computational process extremely 
susceptible to even small numerical errors. 

Near exit plane values of Cy were determined as shown in Figure 14. 
It should be noted that the line of Cy = 0.09, the standard k-£ model 
value, occurring at 0=150° was an assigned default value selected when- 
ever the cross stream velocity gradient became small. Small values of 
Uq r led to extreme variation of the F(/cj^) polynomial, hence the default 
selection in these regions. Pourahmadi and Humphrey [9,10] also note the 
ill behavior of Cy in the core region of the d uct. They attribute this 
behavior to a breakdown in the assumption that u^uj/k is constant in this 
region. Their reasoning is less than compelling. It is entirely possible 
that the problem of determining reasonable values for Cy is completely 
numerical, and that consideration of a more primitive implicit form of 
the Cy relation might yield stable results across the entire duct cross 
section. 

Recognizing that the analysis 3 results are only intermediate, 
lacking a completely acceptable method of selecting C^ , it is useful and 


XXXV-18 







XXXV-20 



interesting to observe position dependent characteristics of the vari- 
ant Cjj . In Figure 14, for 0 = 174° , Cy is observed to take on values 
substantially less than the standard value of 0.09 m near wall regions, 
rising to nearly 0.09 just outside the convex wall region and then fall- 
ing to near zero in the flow core, rising again to near the standard 
value and well beyond just outside the concave wall region. Figure 15 
displays the near concave wall behavior of Cy throughout the flow domain. 
It is evident that the near wall node value of Cy remains essentially 
constant at around 0.03 throughout while the adjacent node value of Cy 
took on default values near the TAD entrance, rising dramatically to a 
value of nearly 0.4, and then falling steadily to a value in excess of 
0.1 before rising dramatically near the duct exit. The turbulent kinet- 
ic energy profiles displayed in Figure 16 exhibit the same qualitative 
character as the unmodified k-£ model profiles of Figure 6 except that 
near outer wall k values are greatly magnified and exhibit a knee as the 
wall is approached. Details of the near convex wall comparisons of tur- 
bulent kinetic energy are shown in Figure 17. Both the standard and 
modified k-£ models predict the same k values at the adjacent wall node, 
however, modified model kinetic energy predictions were somewhat smaller 
than standard at the neighboring y node near the duct exit. Estimated 
velocity and pressure profiles obtained using the modified k-£ model 
were virtually identical to standard model results. 

Two separate BFC computational analyses were performed in order 
to explore the elliptic character of the TAD flow field. These are 
described in Table I, analyses 4 and 5. Before describing results of 
these studies, several general observations are warranted. First, the 
early version of the SSME TAD geometry depicted in Figure 1 suggests a 
computational model with a short entrance region and a relatively long 
straight exit region. In Figure 3 this would indicate flow entrance at 
plane EE and exit at plane AA. Several attempts at such an analysis 
using the BFC option of the PH0ENICS computer code were unsuccessful. 

This does not necessarily indicate a shortcoming in the PH0ENICS/BFC 
solution technique. It is likely that more extensive experience with 
judicious selection of the relaxation factors and other heuristic tuning 
parameters inherent m the PH0ENICS code, and indeed any CFD code with 
like capabilities, would have led to solution convergence with the flow 
orientation described above. Because of the limited time available to 
gain computational experience, and because it was found that extension 
of the straight duct inlet zone greatly enhanced convergence character- 
istics, the suggested SSME TAD flow orientation was reversed, with the 
inlet prescribed at plane AA and exit at plane EE. 

A second observation concerning the BFC computations involves 
the degree of solution convergence at analysis termination. Using the 
absolute residuals sum criterion described more fully in the PH0ENICS 
documentation, the convergence rate at termination had slowed to approx- 
imately 0.1% of total absolute residuals sum per solution sweep. This 
rate was achieved after approximately 30 clock hours of computation on 
the Perkin-Elmer 3220 computer system utilized during this study. 
Although clock run time is a poor absolute measure of computational 
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Figure 16. Turbulent kinetic energy profiles for alaple TAD geometry using the modified k-« 
model: Analysis 3- 
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Figure 17. Comparison of near wall turbulent kinetic energy predictions. 
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efficiency, relative to simple TAD configurations, the BFC analyses re- 
quired an order of magnitude more clock run time and at best achieved 
an absolute residuals sum convergence measure a full order of magnitude 
larger than the simple TAD analyses. Noting that the BFC grid patterns 
were approximately three times the size of the simple cylindrical system 
TAD grids, it is obvious that the BFC analyses were relatively high 
overhead computations. In addition, because calculation was interrupted 
before achieving a proportionately small convergence measure, results of 
the BFC analyses are considered to yield only qualitative predictions of 
flow field characteristics. 

Figures 18, 19, and 20 display results of the analysis 4 BFC com- 
putations for TAD inlet/outlet velocity, pressure, and turbulent kinetic 
energy profiles respectively. These calculations were performed using 
the full grid pattern displayed in Figure 3 with flow inlet at plane AA. 
It is evident in Figures 18 and 19 that BFC calculated TAD inlet flow 
field features are radically different from the inlet boundary assump- 
tions of the simple duct analyses. The inlet velocity profile resembles 
the 1/r proportionality of ideal fluid flow except in the near wall re- 
gions while the outlet profile appears very much like a fully developed 
turbulent duct flow. Obviously the antisymmetric character of the pre- 
dicted TAD outlet pressure profile is quite different from the assumed 
uniform exit plane distribution of the simple geometry analyses. Be- 
cause of the centrifugally induced cross-stream pressure gradient at the 
TAD exit, the duct like flow profile at plane CC developed into the 
stretched, r proportional profile reminiscent of rigid body rotation in 
the 2xD straight exit duct. Predicted developing velocity profiles in 
the exit duct are displayed in Figure 21. In the exit section, the pres- 
sure adjusts to the assumed uniform condition at exit plane EE. The 
core of the TAD exit plane k profile shown in Figure 20 resembles the 
core suppressed profile of the simple geometry analysis with plug flow 
inlet conditions displayed in Figure 6. No inner wall knee is observed 
in the BFC k predictions, however, outer wall turbulent kinetic energies 
are reduced in the BFC analysis. 

In order to reinforce the extended geometry predictions of analy- 
sis 4, a comparison test was devised. BFC analysis 5 was run on the re- 
duced configuration between planes AA and CC. Plug flow inlet conditions 
were prescribed at plane AA and uniform exit plane pressure was assumed. 
Results of this computational procedure at the TAD entrance plane BB were 
used as boundary conditions for a simple geometry analysis using a cylin- 
drical system with the same uniform exit plane pressure assumed. Pre- 
dictions of analysis 5 and 6 should have been identical, at least theo- 
retically. Unfortunately, there was a substantial discrepancy in the 
predicted flow field characteristics as demonstrated in Figures 22 and 
23. Even considering the qualitative nature of the BFC analysis, the 
comparisons were disturbing. Cylindrical system predictions of cross- 
stream velocity gradients at the exit plane were generally much larger 
than results using BFC analysis as is evident in Figure 22. In addi- 
tion, a rather large recirculation zone was predicted using the tradi- 
tional cylindrical coordinate system analysis. Inner wall flow 
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separation was predicted approximately 24° from the exit plane, with the 
counter flow velocity zone extending into the fluid over 7.5% of the duct 
width at the exit plane. By contrast, the near inner wall node velocity 
at the exit plane, obtained using BFC conputations was over 40% of the 
average velocity. Pressure profiles from the two analyses at the TAD 
entrance plane were not even similar m character, as displayed in Fig- 
ure 23. 


If the standard cylindrical system analysis is presumed correct, 
based on an extensive history of computational success and better con- 
vergence characteristics, then the reliability of the current PHOENICS 
BFC option must be suspect. This conclusion must, however, be tempered 
by the fact that acceptable convergence was attained with cylindrical 
system analysis while only run time constrained convergence was achieved 
using BFC’s. In addition, in order to prevent divergence of the 
PHOENICS based computational procedure involving BFC’s, it was necessary 
to neglect certain curvature contributions in the finite difference ap- 
proximations. The rationale for neglect of these terms is more fully 
described in the PHOENICS BFC documentation [21], however, for curvature 
dominated flows, such as that occurring m the SSME TAD, neglect of such 
contributions is dissatisfying. Only by further computational study and 
experimental comparison can deficiencies in the approximation procedure 
be adequately quantified. 
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CONCLUSIONS 


1. The curvature modified k-e turbulence model developed by Pourahmadi 
and Humphrey [9,10] cannot be utilized as presented due to algebraic 
errors in the F(/C^) rooting scheme. 

2. A corrected version of the Pourahmadi and Humphrey [9,10] curvature 
modified turbulence model has potential if an adequate method of 
selecting the multivalued root can be determined. 

3. Neglect of curvature related terms within the BFC option of the gen- 
eral purpose CFD code PHOENICS yields results inconsistent with the 
predictions of standard orthogonal system analyses for sharply 
curved internal flows. 

4. Construction of an acceptable turbulence model for use in the anal- 
ysis of sharply curved internal flows must await the availability of 
adequate experimental data, and will require a systematic program 

of computational verification. 
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